Kinetics of the Melting Transition in DNA 
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We investigate the kinetics of the DNA melting transition using modified versions of the Peyrard- 
Dauxois-Bishop and Poland- Scheraga models that include long and short range interactions. Using 
Brownian dynamics and Monte Carlo simulations, we observe metastable states prior to nucleation 
and demonstrate that the profile and growth modes of the critical droplet can have both classical 
and spinodal characteristics depending on the interaction range and the temperature quench depth. 

PACS numbers: 



The melting of DNA from double to single-stranded 
form has been studied extensively for over forty years Q 
2]. This transition continues to generate interest be- 
cause it provides insight into the biological mechanisms 
of transcription and replication. Although much theo- 
retical work has centered on the role of large nonlinear 
excitations asaprecursor to noeitioe BQfl, surprise 
little focus has been given to the kinetics of this peculiar 
one-dimensional phase transition. 

Nucleation, the process by which a system decays from 
a metastable state, plays an important role in many sys- 
tems undergoing a phase transition During nucle- 
ation, a critical droplet overcomes a free energy barrier 
and grows, causing the system to decay into the stable 
state Q ■ DNA melting is believed to undergo such a pro- 
cess because its sharp melting curve indicates a first-order 
phase transition Q] with the possibility of metastable 
states and because possible nucleation bubbles have been 
observed via electron microscopy 

n 

Although the effect of including helici ty lOfl , nonlinear- 
ities 0,Q], sequence , and defects Il2| in models of 

DNA has been studied extensively, little has been done 
to study how long-range interactions affect these mod- 
els. It is well known that systems with long-range in- 
teractions undergoing phase transitions can be quenched 
into metastable states near a pseudospinodal U. ,, 17|. 
In the mean- field limit, this pseudospinodal becomes a 



well-defined spinodal Q, 0] , which is the limit of 
metastability. Systems undergoing spinodal nucleation 
are driven to the stable phase by diffuse fractal-like 
critical droplets whose amplitudes differ little from the 
metastable background (lslllfll . In contrast, classical 
nucleation is initiated by compact droplets that resemble 
the stable phase. Hence, the inclusion of long-range in- 
teractions in models of DNA has implications for the na- 
ture of nucleation. It is also known that nearest-neighbor 
models of DNA typically underestimate the opening 
probability for small loops 

HQ 

, and experimental ev- 
idence suggests the interactions may be long-range [2l| . 
In this Letter, we analyze how the range of interactions 
affects the kinetics and nature of the DNA melting transi- 
tion in modified versions of the Peyrard-Dauxois-Bishop 
(PDB) 0,0 and Poland-Scheraga (PS) 2] models. 

In our modified version of the PDB model, the state 
of each base pair (bp) is specified by its separation y n 
and velocity y n . The Hamiltonian is given by -ffpDB = 
YaLi m Vi /2 + V(yi)i where the first term is the kinetic 
energy of bps with mass m — 300 amu and the second 
term is a potential given by 

i-l 

V(yi) = D(e~ ayi — l) 2 + W(y i>yj ) 

j=i-R (1) 

- h(y l - y c )(T - T c )Q(yi - y c ). 

The first term in Eq. (JJJ is an on-site Morse potential 
with dissociation energy D = 0.04 eV and a parameter 
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a = 4.45 A" 1 that goes as the inverse well width. The 
Morse potential represents the attraction due to hydro- 
gen bonds and the repulsion of negatively charged phos- 
phate groups. The term W(y n ,y m ) represents the stack- 
ing interaction and the elasticity of the phosphate back- 
bone and is given by the anharmonic potential 

W(y n , y m ) = \K{\ + pe-«(»-+V™)) (y n - y m f. (2) 

Here, the strength of the backbone is given by K = 
0.06 eV/A 2 , and the parameters a = 0.35 A -1 and p = 
0.5 represent the stacking interaction between bound bps. 
The values of these parameters are similar to those in 
the literature In contrast to the PDB model, we 

let our interaction extend beyond nearest-neighbor to all 
bps within a range R, which is allowed to extend to sizes 
of order the persistence length of double-stranded DNA. 
The persistence length, defined as the scale over which 
correlations in the tangential direction of the chain are 
lost, is of the order 100-200 bps in double-stranded DNA. 

The last term in Eq. (JTJ has been added to correct for 
cntropic effects. Simulations of the PDB model do not 
see separation of the strands in a reasonable time because 
the dissociated chain acts as a pseudorandom walker in 
one dimension, which without including higher spatial di- 
mensions quickly returns to the origin even when running 
at temperatures well above the melting temperature. Be- 
cause this effect is not physical, the last term has been 
added to mimic the decrease in the free energy at large 
separations that would be felt by real DNA melting in 
solution. The step function 0(j/j — y c ) turns on the po- 
tential when bps reach separations greater than a critical 
value y c — 2.5 A, above which the bases are assumed to 
be in the stable open phase. The interaction strength 
is given by the parameter h — 0.01 eV/A and by the 
difference between the temperature T and the assumed 
melting temperature T c = 350 K. 

In the PS model, different statistical weights are given 
to bound and unbound segments. A bound segment is 
energetically favored because of hydrogen bonding and 



stacking interactions. Unbound segments are entrop- 
ically favored because the single-stranded segments of 
a loop have a much shorter persistence length allowing 
them to sample a larger configuration of phase space. 
Our modified PS model is described by the Hamiltonian 

N , i 

-§£(2±i) £ (2±i) (3, 

i—1 j—i—R 

- T 5>(4)> 

loops 

where Oi = +1 and Oi — — 1 represent bound and open 
bps, respectively. Here, the binding energy Eq is assumed 
to be the same for all bps and interaction between bps 
Ko extends over a range R. In order to have roughly 
equivalent parameters in both the PDB and PS models, 
we have set Eq = D and Kq — K(l + p). 

The final term in Eq. J2J represents the effective poten- 
tial due to the entropic effects of loops of size I. Using a 
nearest-neighbor droplet model and the above values for 
Eq and Ko, it can be shown that choosing s — 74.4 leads 
to a temperature of 350 K, consistent with our choice 
in the PDB model. In nearest-neighbor PS models, the 
exponent c determines the order of the transition. We 
have chosen c = 2.15, consistent with results for self- 
avoiding loops published previously 13]. The cooper a- 
tivity SI = 0.3 is larger than values reported elsewhere 
because we wish to consider small loops. 

We simulate the PDB and PS models with Brownian 
dynamics (BD) and the Metropolis Monte Carlo (MC) 
algorithm, respectively. In BD the system evolves via a 
Langevin equation 

fc = -W(|fc)-7i/ + !7(t) (4) 

where the noise rj(t) is random Gaussian with <rj(t)> = 
and <7](t)rj(t')> = 2-fk B TS(t - t'). Using a time unit 
t = 1.018 x 10~ 14 s, we choose a time step of 0.25 r and 
a damping constant 7 = 0.1t . In the MC Simula- 
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FIG. 1: Existence of metastable states. The time evolution 
is shown for <y> in the PDB model for (a) PS with R = 1 
and T = 360, (b) PDB with R = 1 and T = 365, (c) PS 
with R = 205 and T = 517, and (d) PDB with R = 205 and 
T = 685. 



tions, random spins are flipped and the change in energy 
between the original and final states is calculated. Nega- 
tive changes in energy are always accepted, and positive 
changes are accepted with probability exp(— AE/ksT). 
For convenience, both simulations use given periodic 
boundary conditions. 

In Figure ^ we plot the evolution of both runs for long 
and short range interactions. In the PDB model, we 
monitor the growth of the mean separation of the strands 
<y> — jfYl Hi, and for the PS model we observe the 
magnetization M = ■hY^^i' Each run is instantaneously 
raised from a low temperature to a temperature above 
the melting temperature. In PDB we see that for both 
R = 1 and R = 205, the separation stabilizes after the 
quench before growing rapidly, indicating that the system 
enters a metastable state before nucleating. Likewise, 
the PS model exhibits similar metastability in which the 
magnetization M = jj^Pi stabilizes before the system 
nucleates into the stable melted state M = 1. 

In both the PDB and PS models metastability is ob- 
served for all values of R investigated as long as the 
temperature is not too high. In general, systems with 
longer ranges have longer lifetimes and exhibit metasta- 



bility at higher temperatures This behavior is due 
to the fact that the metastable lifetime is proportional 
to the free energy cost of a critical droplet, which in 
a one-dimensional system scales linearly with the range 
Long-range systems are closer to mean- 
field and can be quenched to temperatures close to the 
mean-field spinodal, while short range systems reach the 
Becker-Doring limit where the system nucleates before 
metastable equilibrium can be achieved. We have cho- 
sen temperatures for each range such that metastability 
can be observed in a reasonable time. In both models 
we find that the nucleation rate decreases with increas- 
ing range and that systems with longer ranges can reach 
metastable equilibrium at higher temperatures than at 
shorter ranges. 

Nucleation can be viewed as passing from a metastable 
well to a stable well in a free ^landscape Q. The 
critical droplet is the system configuration at the top of 
a saddle point ridge where the probability of decaying 
into the metastable state equals the probability of grow- 
ing into the stable state. In order to determine when this 
droplet occurs during a run, we make multiple copies of 
the system in its initial state and rerun it while adding a 
random perturbation at some intervention time tint [20j. 
If after the perturbation the system is more likely to de- 
cay to the metastable state, then we have not yet reached 
the critical droplet, while if it tends to nucleate than the 
droplet is past critical and has entered the growth phase. 
When the probability of returning to the metastable state 
and growing to the stable phase are equal, we say the con- 
figuration is the critical droplet. In the PDB model, we 
intervene by turning off the noise at tint and observing 
whether the system grows or decays. In the MC simu- 
lations of the PS model, we change the random number 
sequence at tint and rerun the system multiple times to 
see when the probability of nucleation is roughly one half. 

The profiles of the critical droplets are shown in Fig- 
ure El To reduce the noise, each profile has been coarse- 
grained by averaging over all bps within the interaction 
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droplet profiles at various times after the critical droplet 
and subtracting the critical droplet profile. This averag- 
ing is done to reduce the noise in the growth mode. For 
both models, small R and T gives growth modes that 
are largest at the surface of the droplet, and large R and 
T systems grow mostly in the center. This difference is 
expected as classical droplets resemble the stable phase 
and grow at the surface, and diffuse spinodal droplets 
initially grow by filling in at the center [19(. 



FIG. 2: (color online). Plots of the critical droplet (blue 
squares), growth mode (red triangles), and the fit (green line) 
for (a) PS with R = 1 and T = 360, (b) PDB with R = 1 and 
T = 365, (c) PS with R = 205 and T = 517, and (d) PDB 
with R = 205 and T = 685. 

range. In both models, smaller R and T give com- 
pact droplets as expected in classical nucleation. Runs 
with R » 1 and high T give diffuse, small amplitude 
droplet profiles that fit well to an inverse hyperbolic co- 
sine squared shape. These results are consistent with 
spinodal nucleation [is| . 

We also plot in Figure |3 the growth modes of the 
droplet. Growth modes are obtained by averaging the 



In summary, we have observed that nucleation drives 
the melting transition in modified PDB and PS mod- 
els of DNA. We find that long-range interactions, for 
which there is experimental evidence, and deep quench 
depths give diffuse droplets reminiscent of spinodal nu- 
cleation while the short-range shallow quenches give com- 
pact classical droplets. These results in addition to ex- 
perimental evidence for long-range interactions [2 1| in 
real DNA suggest the possibility that nucleation in DNA 
may have spinodal characteristics. 
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